import numpy as np
import matplotlib.pyplot as plt
import matplotlib.pyplot as plt1
import matplotlib.pyplot as plt2
import matplotlib.pyplot as plt3
import matplotlib.image as mpimg
from PIL import Image  # 需安装：pip install pillow
from matplotlib import cm
from scipy.interpolate import griddata
import matplotlib.colors as colors
from mpl_toolkits.axes_grid1 import make_axes_locatable
from matplotlib.patches import Rectangle
import io
import json
from scipy import signal
import pickle
import os
from io import BytesIO  # 内置模块，无需额外安装
from decimal import Decimal

temp = []
IPValue = {}
# fig, ax = plt.subplots(figsize=(15, 7.5))
# plt.figure(figsize=(15, 7.5))
# 倒置y轴
#plt.gca().invert_yaxis()
# fig.set_size_inches(15, 7.5)  # 设置大小为10x5英寸
# 全局设置默认画布大小
#plt.rcParams['figure.figsize'] = [15, 7.5]  # 所有图形默认大小为10x5英寸
#fig1, ax1 = plt1.subplots(figsize=(15, 7.5))

def ClosePlt():
    plt.close()

def ERTuiArray(ABMN, a1_dtype, a1_shape, Aspacing):
    abmn = ABMN.encode('ISO-8859-1')
    ABMN1 = np.fromstring(abmn, dtype=getattr(np, a1_dtype)).reshape(eval(a1_shape))
    # print(print(ABMN1[0,:]))
    return ABMN1

#拟断面开始
def SectionArray(jb_tu,json_str, ns,nss,types,
                 obj1,freq_Indx, rx_voltage_data_len,num_points,num_points1,filter,eleNum,model,CI):
    # 将JSON字符串转换为Python对象
    data = json.loads(json_str)
    # ReadData(data)
    io = ERTuiSection(data,jb_tu, int(ns), ns,nss,types,obj1,freq_Indx, rx_voltage_data_len,num_points,num_points1,filter,eleNum,model,CI)
    return io


def ERTuiSectionNs(json_str):
    data = json.loads(json_str)
    # print("json abmn =======================================")
    # print(data)
    # NS = np.sum(range(1,int(data['Max_n_spacing'])+1)) + int(data['Max_n_spacing'])*((int(data['num_electrodes'])-3)-int(data['Max_n_spacing']))      #总排列数
    NS = len(ERTuiArray(data['ABMN'], data['dtype'], data['shape'], data['Aspacing']))
   # print(NS)
    return int(NS)

def ERTuiSection(data,jb_tu, AI, ns,nss,types,obj1,freq_Indx, rx_voltage_data_len,num_points,num_points1,filter,eleNum,model,dtypelen0):
    # NS = np.sum(range(1,int(data['Max_n_spacing'])+1)) + int(data['Max_n_spacing'])*((int(data['num_electrodes'])-3)-int(data['Max_n_spacing']))      #总排列数
    ###3. 读取CVS数据及系统参数
    '''
        手部首先需要向主机发送数据采集的命令，规定采集的参数，然后主机才能返回数据。这里假设命令发送成果，数据成功上传。
        数据按每次叠加上传，例如我们设置命令为叠加5次，每次波形长度为8，那么手部会收到8个完整波形后再做相应的处理。
        这8个波形的数据在数据文件中是一个数据块，有相应的头文件与之对应。
    '''
    Data = jb_tu
    ABMN = ERTuiArray(data['ABMN'], data['dtype'], data['shape'], data['Aspacing'])
    # list1 = list(map(int,data['ABMN']))
    #
    # ABMN = np.array(list1)
    # IPValue = {}
    global IPValue
    IPValue = IPValue;
    #plt.ion()  # 打开交互模式
    # for AI in range(NS):
    #Data = ReadData(data, ABMN[AI, :])
    #print("----------------------ABMN[AI,:]------------")
    # print(ABMN)
    #print(AI)
    data['ABMN'] = ABMN
    ###4. 处理数据

    IPValue[AI] = ComputeRhoIP(Data, data, AI,rx_voltage_data_len)
    #print(IPValue)
    ###5. 画出图形
    io = []
    io1 = ERTPlotSection(IPValue, ABMN, int(data['Aspacing']), ns,nss,eleNum,model,dtypelen0)
    #io2 = ERTPlotArray(IPValue, ABMN, int(data['Aspacing']))
    # io3 = draw_fft_square_wave(obj1, freq_Indx, rx_voltage_data_len)
    # io4 = draw_fft(num_points, num_points1)
    NIP = len(IPValue)

    if types == 0: #拟断面
        #plt.close()
        io = []
        #io = io1
        io.append(io1)
    elif types == 1: #排列图
        #plt1.close()
        io = []
        iop = ERTPlotArray(IPValue, ABMN, int(data['Aspacing']),eleNum,model,dtypelen0)
        io.append(iop)
        io.append(io1)
    elif types == 2: #方波
        #plt.close()
        io = []
        iof = draw_fft_square_wave(obj1, freq_Indx, rx_voltage_data_len,filter)
        io.append(iof)
        io.append(io1)
    elif types == 3: #频谱
        #plt1.close()
        #io = draw_fft(num_points, num_points1)
        io = []
        iopin = draw_fft(obj1, freq_Indx, rx_voltage_data_len,filter)
        io.append(iopin)
        io.append(io1)
        print("types == 3: #频谱: ",len(io))
    if nss == NIP:
        io = []
        io.append(io1)

        io.append(ERTPlotArray(IPValue, ABMN, int(data['Aspacing']),eleNum,model,dtypelen0))
        io3 = draw_fft_square_wave(obj1, freq_Indx, rx_voltage_data_len,filter)
        #io4 = draw_fft(num_points, num_points1)
        io4 = draw_fft(obj1, freq_Indx, rx_voltage_data_len,filter)
        io.append(io3)
        io.append(io4)
        plt.close();
        IPValue = {}
        return io
    else:
        ios=[]
        #print("===============================================================================")
        ios.append(io)
    return io

def ComputeRhoIP(Data, data, AI, NFFT):
    '''
        计算电阻，极化率，及画图的位置.
        此算法暂时只针对偶极-偶极排列
        :param Data:
        :param SysInfo:
        :param SurveyInfo:
        :param ABMN:
        :return: IPValues
    '''
    print("-------------------------------------------"+str(NFFT))
    # ABMN1 = ERTuiArray(data['ABMN'], data['dtype'], data['shape'], data['Aspacing'])
    # ABMN = ABMN1[AI, :]
    ABMN1 = data['ABMN']
    ABMN = ABMN1[AI, :]
    # print("----------------------------------------ABMN[AI,:]------------")
    # print(ABMN)
    SR = int(data['RXSR'])  # 数据点
    TXFrq = 0
    Mpts = 0
    Mst =0
    VpPnts =0;
    if data['TXFrq'] == '1/2':
        TXFrq = 0.5  # 发射频率
        Mpts = 167
        Mst = 118
        VpPnts= 128
    elif data['TXFrq'] == '1/4':
        TXFrq = 0.25  # 发射频率
        Mpts = 332
        Mst = 234
        VpPnts= 256
    elif data['TXFrq'] == '1/8':
        TXFrq = 0.125  # 发射频率
        Mpts = 666
        Mst = 464
        VpPnts= 512
    else:
        TXFrq = int(data['TXFrq'])  # 发射频率
        if data['TXFrq'] == '1':
            Mpts = 83
            Mst = 61
            VpPnts= 64
        elif data['TXFrq'] == '2':
            Mpts = 1331
            Mst = 925
            VpPnts= 1024
        elif data['TXFrq'] == '4':
            Mpts = 666
            Mst = 464
            VpPnts= 512
        elif data['TXFrq'] == '8':
            Mpts = 332
            Mst = 234
            VpPnts= 256

    print("++++++++++++Mpts Mst VpPnts+++++++++++++++++")
    print(Mpts,Mst,VpPnts)
    StackNO = int(data['StackNO'])  # 叠加次数  这三个都是app填写的
    Aspacing = int(data['Aspacing'])
    TXA = (int(float(ABMN[0])) - 1) * Aspacing  # Aspacing 点距 也是app 填写的
    TXB = (int(float(ABMN[1])) - 1) * Aspacing
    RXM = (int(float(ABMN[2])) - 1) * Aspacing
    RXN = (int(float(ABMN[3])) - 1) * Aspacing


    IPValues = {}
    IPValues['TXA'] = TXA
    IPValues['TXB'] = TXB
    IPValues['RXM'] = RXM
    IPValues['RXN'] = RXN
    IPValues['X'] = (RXM - TXB) / 2 + TXB  # 画图横坐标位置
    IPValues['Z'] = (RXM - TXB) / 2 + Aspacing / 2  # 画图垂向坐标位置
    IPValues['n'] = ABMN[2] - ABMN[1]  # 隔离系数

    # print("--------------IPValues['n']---------------")
    # print(IPValues['n'])

    #print(type(Data))
    #datas = Data.replace("'"," ");
    #dict = json.loads(s=Data)
    dict = eval(Data)
    # 使用字典推导式将字符串键转换为整数键
    Data = {int(key): value for key, value in dict.items()}
    #print(type(Data),Data)
    # 将字典中的值转换为浮点数类型
    for key in Data:
        # new_dict[key] = float(eval(new_dict[key]))
        Data[key] = np.array(Data[key])
   # print("-----------Data --------")
    #print(type(Data),Data)
    # 叠加处理数据
    temp =  Data[0]  # 将第一个叠加数据暂存起来
    # print(type(temp),temp)
    for BI in range(1, StackNO):
        #print(StackNO,BI)
        temp = temp + np.array(Data[BI])  # 叠加数据
    temp = temp / StackNO  # 平均数据
    #print(type(temp),temp)
    #print("---------------- SR ="+str(SR))
   #print("---------------- TXFrq ="+str(TXFrq))
    #NP = int(SR/TXFrq)  # 一个周期的数据点数
    NP = SR  # 一个周期的数据点数
    #print("---------------- np =")
   # print(NP)
    NPhlf = int(NP / 2)  # 半个周期数据点数
    #print("---------------- NPhlf ="+str(NPhlf))
    tmp = (temp[0:NPhlf, :] - temp[NPhlf:, :]) / 2  # 将后半个周期和前半个周期平均
    # print("---------------- tmp ="+str(tmp))
    TX_I = tmp[:, 0]  # 发射电流数据
    # print("---------------- TX_I ="+str(TX_I))
    RX_V = tmp[:, 1]  # 接收电压数据
    #VpPnts = int(0.0625 * NP)  # 数据点数用于计算电流和电压
    #TXI = np.sum(TX_I[int(NPhlf / 2 - VpPnts):int(NPhlf / 2)]) / VpPnts  # 发射电流值
    # Vp = np.sum(RX_V[int(NPhlf / 2 - VpPnts):int(NPhlf / 2)]) / VpPnts  # 接收电压值
    #FFT转化
    FFT =  2 * np.fft.fft(TX_I, n=NFFT)/NFFT
    NFFT_calc = int(NFFT/2) + 1
    TXI =  np.sum(np.abs(FFT[0:NFFT_calc]))

    FFT_R =  2 * np.fft.fft(RX_V, n=NFFT)/NFFT
    NFFT_R_calc = int(NFFT/2) + 1
    Vp =  np.sum(np.abs(FFT_R[0:NFFT_R_calc]))

    IPValues['TXI'] = TXI
    IPValues['RXV'] = Vp
    # Mpts = int(0.08125 * NP)
    # Mst = int(0.05625 * NP) + 4
    # print()
    IP = 1 / TXFrq / NP * 1.87 / Vp * np.sum(RX_V[int(NPhlf + Mst):int(NPhlf + Mst + Mpts)]) / 10
    IPValues['IP'] = IP
    #print("----------Rho-----------")
    print(Vp,TXI,np.pi,IPValues['n'],Aspacing)
    #print(Vp,TXI)
    #Rho = 2.0/2.0 * np.pi * IPValues['n'] * (IPValues['n'] + 1) * (IPValues['n'] + 2) * Aspacing
    Rho = Vp/TXI * np.pi * IPValues['n'] * (IPValues['n'] + 1) * (IPValues['n'] + 2) * Aspacing
    print("------------ Rho ------------"+str(Rho))
    IPValues['Rho'] = Rho

    #print(">>>>>>>>>>>>>>>>>>>>>IPValues= ")
   # print(IPValues)
    return IPValues


def ERTPlotSection(IPValues, ABMN, Aspacing, NS,NSS,eleNum,model,CI):
    #print("___________________ERTPlot__________________11111")
    '''
        本程序画出高密度电阻率法的等值线图和排列图
        :param IPValues:
        :return: NONE
    '''
    ####这里画出两个图，一个是排列图，利用红色圆圈表示发射，用黑色圆表示接收
    # 最后一个IP的值
    # 获取相应的画图参数
    # 设置画布大小为 10x6 英寸
    # 1. 定义沙盒图片路径（替换为实际包名和图片名）
    package_name = "com.example.shiwu"  # 你的应用包名
    image_name = "1-1c6bf17a313a4a5a806454fd7d9f0fec_niduan.png"  # 沙盒中的图片文件名
    sandbox_image_path = f"/data/data/{package_name}/files/img/{image_name}"
    #img = mpimg.imread(sandbox_image_path)  # 读取图片
    #plt.imshow(img)
    try:
        img = Image.open(sandbox_image_path)
    except FileNotFoundError:
        print("错误：文件不存在")
        fig, ax = plt.subplots(figsize=(15, 7.5))
        #ax.imshow(img_array)  # 显示原图作为底图
        plt.figure(figsize=(15, 7.5))
        print("-------------------------IPValues---------------------------")
        print(type(IPValues))
        NIP = len(IPValues)
        print("-------------------------NIP---------------------------")
        print(NIP)
        if NIP == 1:
            lastIP = IPValues[NIP-1]
        else:
            last_key = next(reversed(IPValues))
            lastIP = IPValues[last_key]

        X = lastIP['X']
        Z = lastIP['Z']
        rho = lastIP['Rho']
        # #####计算图的尺寸
        #print(ABMN[:, 0])
        Min_X = np.min(ABMN[:, 0] - 1) * Aspacing
        # Max_X = np.max(ABMN[:, 3] - 1)+2
        # Max_X = np.max(ABMN[:, 3] - 1) * Aspacing
        # Min_Z = np.min(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
        # Max_Z = np.max(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
        if model == 0 or model == 1:
            Max_X = np.max(ABMN[:, 3])+2
            Min_Z = int(eleNum)
            Max_Z = int(eleNum)*Aspacing
        else:
            Max_X = np.max(ABMN[:, 3])+2
            Min_Z = int(eleNum)
            Max_Z = int(eleNum)*Aspacing
        BE_MAX = np.max(ABMN[:,3])
        #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
        A_indx_sum = np.sum(ABMN[:,0] == 1)+1
        #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
        A_indx_sum_ouji = np.sum(ABMN[:,3] == BE_MAX)+1

        # fig1 = plt.figure(num=1)
        #
        # ax = fig1.add_subplot(111)
        # ax.legend(frameon=False)

        colormap = cm.get_cmap('jet')
        rho_range = np.logspace(1, 5, 40)

        # plt.xlim([Min_X, Max_X])
        # plt.ylim([0, Max_Z + 10])
        plt.xticks(np.arange(Min_X, Max_X, 1))
        if model == 0 or model == 1:
            abs_yticks=np.arange( -(A_indx_sum_ouji-1),0)
            print("abs_yticks = 0/1")
            print(abs_yticks)
            print(abs_yticks)
            plt.yticks(abs_yticks)
        else:
            abs_yticks=np.arange(-(A_indx_sum-1),0 )
            print("abs_yticks = 3")
            print(abs_yticks)
            plt.yticks(abs_yticks)
            # plt1.yticks(np.insert(np.arange(Min_Z, Max_Z + Aspacing / 2, Aspacing / 2), 0, 0))

        NS = np.size(ABMN, 0)

        ###画出所有的IPvalues
        rho_trans = (np.log10(rho) - np.log10(min(rho_range))) / 4
        cls = colormap(rho_trans)
        # 画为圆或者方块都可以，这里选择了方块，使得两个点之间的电阻率图像可以连接起来
        # Rho_C = plt.Circle((X, Z), radius=Aspacing / 3, fill=True, color=cls, linewidth=2.0, \
        #                   clip_on=False)   # 拟断面图画图点
        # plt.gca().add_artist(Rho_C)
        x1=1
        y1=1
        if model == 0 or model == 1:
            BE = int(ABMN[CI, 1])
            ME = int(ABMN[CI, 2])
            XE = int(ME - BE) / 2 + BE
            ZE = -int(ME - BE)

            Rho_rect = Rectangle((XE-x1/2, ZE-y1/2), x1, y1,
                                 fill=True, color=cls, \
                                 clip_on=False)  # 拟断面图画图点
        else:
            ME = ABMN[CI, 1]
            NE = ABMN[CI, 2]
            BE = ABMN[CI, 3]
            XE = (NE - ME) / 2 + ME
            ZE = -(BE - NE)
            Rho_rect = Rectangle((XE-x1/2, ZE-y1/2), x1, y1,
                                 fill=True, color=cls, \
                                 clip_on=False)  # 拟断面图画图点
        plt.gca().add_patch(Rho_rect)  # 绘制图像形状方形
        print("XE, ZE ",XE, ZE)
        plt.annotate(format(rho, '.0f'), xy=(XE, ZE), fontsize=8, ha='center', va='center')
        #plt.gca().invert_yaxis()
        plt.gca().set_xlabel('Array Location(m)', size=14)
        plt.gca().set_ylabel('n-spacing', size=14)
        plt.gca().xaxis.set_label_position('top')
        plt.gca().xaxis.tick_top()
        plt.gca().set_aspect('equal')
        # 用隔离系数代替纵向坐标的标签
        labels = [item.get_text() for item in ax.get_yticklabels()]
        # print(labels)
        #NLAB = len(labels)
        # if type == 0 or type == 1:
        #     labs = np.arange(-(A_indx_sum_ouji-1),1)
        # else:
        #     labs = np.arange(0, A_indx_sum, 1)

        # abs_labs = [abs(num) for num in labs]
        # abs_labs = labels[::-1]
        # print(abs_labs)
        ax.set_yticklabels(abs_yticks)
        ax.spines['bottom'].set_visible(False)
        # labels = [item.get_text() for item in ax.get_yticklabels()]
        # NLAB = len(labels)
        # labs = np.insert(np.arange(int(np.min(ABMN[:, 2] - ABMN[:, 1])), \
        #                            int(np.max(ABMN[:, 2] - ABMN[:, 1]) + 1), 1), 0, 0)
        #print("-----------------------------------NLAB--------------------------------------------")
        # print(NIP)
        # print(NSS)
        if NIP == NSS:
            plt.clf()
            Max_X = np.max(ABMN[:, 3] - 1) * Aspacing
            Min_Z = np.min(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
            Max_Z = np.max(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
            ax.axes.get_xaxis().set_visible(False)
            CX = np.zeros((NSS,))
            CZ = np.zeros((NSS,))
            Rho = np.zeros((NSS,))
            for I in range(NSS):
                currentIP = IPValues[I]
                CX[I] = currentIP['X']
                CZ[I] = currentIP['Z']
                Rho[I] = currentIP['Rho']
            # print("CX, CZ , Rho")
            print("CX, CZ , Rho",CX, CZ , Rho)
            #生成网格
            xv, zv = np.meshgrid(np.linspace(Min_X, Max_X, 100), np.linspace(Min_Z, Max_Z, 100))
            #插值
            rho_av = griddata((CX, CZ), Rho, (xv, zv))
            cs = plt.contourf(xv, zv, rho_av, levels=np.logspace(1,5,40), cmap=cm.jet,norm = colors.LogNorm())
            cs.set_clim(vmin=10, vmax=10000)
            # plt.gca().clabel(cs, fmt='%2.1f', colors='k', fontsize=8, inline=True,inline_spacing=0)


            plt.gca().invert_yaxis()
            plt.gca().set_xlabel('Array Location(m)', size=14)
            plt.gca().set_ylabel('n-spacing', size=14)
            plt.gca().xaxis.set_label_position('top')
            plt.gca().xaxis.tick_top()
            plt.gca().set_aspect('equal')
            plt.yticks(np.insert(np.arange(Min_Z, Max_Z + Aspacing / 2, Aspacing / 2), 0, 0))
            plt.xticks(np.arange(Min_X, Max_X, Aspacing))


            # plt.yticks(abs_yticks)
            # plt.xticks(np.arange(Min_X, Max_X, 1))
            labels = [item.get_text() for item in ax.get_yticklabels()]
            NLAB = len(labels)
            labs = np.insert(np.arange(int(np.min(ABMN[:, 2] - ABMN[:, 1])), \
                                       int(np.max(ABMN[:, 2] - ABMN[:, 1]) + 1), 1), 0, 0)
            if NLAB == len(labs):
                ax.set_yticklabels(labs)
            divider = make_axes_locatable(ax)
            cax = divider.append_axes('right', size="2%", pad=18.7)
            #cax = divider.new_vertical('right',size="5%", pad=0.7, pack_start=True)
            sm = plt.cm.ScalarMappable(norm=colors.LogNorm(), cmap=colormap)
            sm.set_clim(vmin=10, vmax=10000)
            cbar = plt.colorbar(sm, cax=cax)
            cbar.set_label('Apparent Resistivity(Ohm-m)', rotation=270, labelpad=20)
            ax.spines['bottom'].set_visible(False)
            plt.show()
            #cbar.remove()
            f = io.BytesIO()
            plt.savefig(f, format="png")
            fio = f.getvalue()
            #print("___________________ERTPlot__________________33333")
            plt.close()
            f.close()
            return fio

        f = io.BytesIO()
        plt.savefig(f, format="png")
        fio = f.getvalue()
        plt.close()
        f.close()
        # 保存图表对象为pkl文件


        #print("___________________ERTPlot__________________22222")
        return fio
    else:
        print("文件内容：", img)  # 只有无异常时才执行
        # fig, axs = image_to_fig_ax(sandbox_image_path)
        # 加载已有的图表对象
        pkl_path = os.path.join("/data/data/com.example.shiwu/files/img/", "matplotlib_figure.pkl")
        with open(pkl_path, 'rb') as f:
            fig = pickle.load(f)

        # 获取当前轴对象
        axs = fig.gca()
        # print("------------------------IPValues---------------------------")
        # print(type(IPValues))
        NIP = len(IPValues)
        # print("-------------------------NIP---------------------------")
        # print(NIP)
        if NIP == 1:
            lastIP = IPValues[NIP-1]
        else:
            last_key = next(reversed(IPValues))
            lastIP = IPValues[last_key]

        X = lastIP['X']
        Z = lastIP['Z']
        rho = lastIP['Rho']
        # #####计算图的尺寸
        Min_X = np.min(ABMN[:, 0] - 1) * Aspacing
        # Max_Z = np.max(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
        if model == 0 or model == 1:
            Max_X = np.max(ABMN[:, 3])+2
            Min_Z = int(eleNum)
            Max_Z = int(eleNum)*Aspacing
        else:
            Max_X = np.max(ABMN[:, 3])+2
            Min_Z = int(eleNum)
            Max_Z = int(eleNum)*Aspacing
        BE_MAX = np.max(ABMN[:,3])
        #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
        A_indx_sum = np.sum(ABMN[:,0] == 1)+1
        #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
        A_indx_sum_ouji = np.sum(ABMN[:,3] == BE_MAX)+1


        colormap = cm.get_cmap('jet')
        rho_range = np.logspace(1, 5, 40)

        plt.xticks(np.arange(Min_X, Max_X, 1))
        if model == 0 or model == 1:
            abs_yticks=np.arange( -(A_indx_sum_ouji-1),0)
            print("abs_yticks = 0/1")
            print(abs_yticks)
            print(abs_yticks)
            plt.yticks(abs_yticks)
        else:
            abs_yticks=np.arange(-(A_indx_sum-1),0 )
            print("abs_yticks = 3")
            print(abs_yticks)
            plt.yticks(abs_yticks)

        NS = np.size(ABMN, 0)

        ###画出所有的IPvalues
        rho_trans = (np.log10(rho) - np.log10(min(rho_range))) / 4
        cls = colormap(rho_trans)
        x1=1
        y1=1
        if model == 0 or model == 1:
            BE = int(ABMN[CI, 1])
            ME = int(ABMN[CI, 2])
            XE = int(ME - BE) / 2 + BE
            ZE = -int(ME - BE)

            Rho_rect = Rectangle((XE-x1/2, ZE-y1/2), x1, y1,
                             fill=True, color=cls, \
                             clip_on=False)  # 拟断面图画图点
        else:
            ME = ABMN[CI, 1]
            NE = ABMN[CI, 2]
            BE = ABMN[CI, 3]
            XE = (NE - ME) / 2 + ME
            ZE = -(BE - NE)
            Rho_rect = Rectangle((XE-x1/2, ZE-y1/2), x1, y1,
                                 fill=True, color=cls, \
                                 clip_on=False)  # 拟断面图画图点
            # 隐藏坐标轴（可选）
        # axs.axis('off')
        axs.add_patch(Rho_rect)  # 绘制图像形状方形
        # axs.add_patch(Rectangle((100, 100), 200, 200,
        #                         fill=True, color=cls, \
        #                         clip_on=False))


        print("XE, ZE ",XE, ZE)
        plt.annotate(format(rho, '.0f'), xy=(XE, ZE), fontsize=8, ha='center', va='center')
        # plt.gca().invert_yaxis()
        plt.gca().set_xlabel('Array Location(m)', size=14)
        plt.gca().set_ylabel('n-spacing', size=14)
        plt.gca().xaxis.set_label_position('top')
        plt.gca().xaxis.tick_top()
        plt.gca().set_aspect('equal')


        # 用隔离系数代替纵向坐标的标签
        # labels = [item.get_text() for item in ax.get_yticklabels()]

        # axs.set_yticklabels(abs_yticks)
        # axs.spines['bottom'].set_visible(False)
        if NIP == NSS:
            plt.clf()
            Max_X = np.max(ABMN[:, 3] - 1) * Aspacing
            Min_Z = np.min(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
            Max_Z = np.max(ABMN[:, 2] - ABMN[:, 1]) * Aspacing / 2 + Aspacing / 2
            axs.axes.get_xaxis().set_visible(False)
            CX = np.zeros((NSS,))
            CZ = np.zeros((NSS,))
            Rho = np.zeros((NSS,))
            for I in range(NSS):
                currentIP = IPValues[I]
                CX[I] = currentIP['X']
                CZ[I] = currentIP['Z']
                Rho[I] = currentIP['Rho']
            # print("CX, CZ , Rho")
            print("CX, CZ , Rho",CX, CZ , Rho)
            #生成网格
            xv, zv = np.meshgrid(np.linspace(Min_X, Max_X, 100), np.linspace(Min_Z, Max_Z, 100))
            #插值
            rho_av = griddata((CX, CZ), Rho, (xv, zv))
            cs = plt.contourf(xv, zv, rho_av, levels=np.logspace(1,5,40), cmap=cm.jet,norm = colors.LogNorm())
            cs.set_clim(vmin=10, vmax=10000)
            # plt.gca().clabel(cs, fmt='%2.1f', colors='k', fontsize=8, inline=True,inline_spacing=0)
            plt.gca().invert_yaxis()
            plt.gca().set_xlabel('Array Location(m)', size=14)
            plt.gca().set_ylabel('n-spacing', size=14)
            plt.gca().xaxis.set_label_position('top')
            plt.gca().xaxis.tick_top()
            plt.gca().set_aspect('equal')
            plt.yticks(np.insert(np.arange(Min_Z, Max_Z + Aspacing / 2, Aspacing / 2), 0, 0))
            plt.xticks(np.arange(Min_X, Max_X, Aspacing))
            # plt.yticks(abs_yticks)
            # plt.xticks(np.arange(Min_X, Max_X, 1))
            labels = [item.get_text() for item in axs.get_yticklabels()]
            NLAB = len(labels)
            labs = np.insert(np.arange(int(np.min(ABMN[:, 2] - ABMN[:, 1])), \
                                       int(np.max(ABMN[:, 2] - ABMN[:, 1]) + 1), 1), 0, 0)
            if NLAB == len(labs):
                axs.set_yticklabels(labs)
            divider = make_axes_locatable(axs)
            cax = divider.append_axes('right', size="2%", pad=18.7)
            #cax = divider.new_vertical('right',size="5%", pad=0.7, pack_start=True)
            sm = plt.cm.ScalarMappable(norm=colors.LogNorm(), cmap=colormap)
            sm.set_clim(vmin=10, vmax=10000)
            cbar = plt.colorbar(sm, cax=cax)
            cbar.set_label('Apparent Resistivity(Ohm-m)', rotation=270, labelpad=20)
            axs.spines['bottom'].set_visible(False)
            plt.show()
            #cbar.remove()
            f = io.BytesIO()
            plt.savefig(f, format="png")
            #print("___________________ERTPlot__________________33333")
            return f.getvalue()

        # 获取Android应用的私有文件目录（沙盒目录）
        # 该目录路径通常为：/data/data/[应用包名]/files/
        # 构建保存路径（在沙盒目录下创建pkl文件）
        with open(pkl_path, 'wb') as f:
            pickle.dump(fig, f)  # 保存整个Figure对象
            print("图表已保存为 matplotlib_figure.pkl")

        f = io.BytesIO()
        plt.tight_layout()
        plt.savefig(f, format="png")
        fio = f.getvalue()
        plt.close()
        f.close()
        #print("___________________ERTPlot__________________22222")
        return fio
#排列图示
def ERTPlotArray(IPValues, ABMN, Aspacing,eleNum,model,CI):
    # print(eleNum)
    # print(model)
    # print(CI)
    '''
        本程序画出高密度电阻率法的等值线图和排列图
        :param IPValues:
        :return: NONE
    '''
    ####这里画出两个图，一个是排列图，利用红色圆圈表示发射，用黑色圆表示接收
    #最后一个IP的值
    #获取相应的画图参数
    NIP = len(IPValues)
    if NIP == 1:
        lastIP = IPValues[NIP-1]
    else:
        last_key = next(reversed(IPValues))
        lastIP = IPValues[last_key]

    # TXA = lastIP['TXA']
    # TXB = lastIP['TXB']
    # RXM = lastIP['RXM']
    # RXN = lastIP['RXN']
    X = lastIP['X']
    Z = lastIP['Z']
    rho = lastIP['Rho']
    # #####计算图的尺寸
    Min_X = np.min(ABMN[:,0]-1) * Aspacing
    # Max_X = np.max(ABMN[:,3]-1) * Aspacing
    # Min_Z = np.min(ABMN[:,2]-ABMN[:,1])*Aspacing/2 + Aspacing/2
    # Max_Z = np.max(ABMN[:,2]-ABMN[:,1])*Aspacing/2 + Aspacing/2
    Max_X = np.max(ABMN[:, 3])+2
    Min_Z = int(eleNum)
    Max_Z = int(eleNum)*Aspacing
    BE_MAX = np.max(ABMN[:,3])
    #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
    A_indx_sum = np.sum(ABMN[:,0] == 1)+1
    #A_indx_sum 表示ABMN在温纳、斯伦贝谢、温纳-斯伦贝谢装置A=1时出现的次数作为深度的y轴的刻度
    A_indx_sum_ouji = np.sum(ABMN[:,3] == BE_MAX)+1
    # fig1 = plt.figure(num=1)
    #
    # ax = fig1.add_subplot(111)
    # ax.legend(frameon=False)

    colormap = cm.get_cmap('jet')
    rho_range = np.logspace(1,5,40)
    fig1,ax1 = plt1.subplots(figsize=(20, 10))
    # plt1.figure(figsize=(15, 7.5))
    plt1.xlim(0)  # 设置横向尺度
    plt.ylim([ -int(eleNum),0])
    #plt.figure(figsize=(22, 10))
    # 设置纵向尺度
    Radius = 0.1 # 圆表示画图点，圆的半径
    # 画出所有排列的地下记录点
    NS = np.size(ABMN, 0)
    print("lastIP: ",lastIP)


    if model == 0 or model == 1:
        BE = int(ABMN[CI, 1])
        ME = int(ABMN[CI, 2])
        XE = int(ME - BE) / 2 + BE
        ZE = -int(ME - BE)
        print("XE, ZE:",XE, ZE)

        E_C = plt.Circle((XE, ZE), radius=Radius, fill=False, color='black', linewidth=1.0, \
                         clip_on=False)  # 画图点位置
        plt1.gca().add_artist(E_C)
        # 画出本次测量的地下记录的及发射接收位置
        TXA_C = plt.Circle((ABMN[CI][0], 0), radius=Radius, fill=False, color='red', linewidth=2.0, \
                           clip_on=False)  # 发射A位置
        plt.gca().add_artist(TXA_C)
        TXB_C = plt.Circle((ABMN[CI][1], 0), radius=Radius, fill=False, color='red', linewidth=2.0, \
                           clip_on=False)  # 发射B位置
        plt.gca().add_artist(TXB_C)
        RXM_C = plt.Circle((ABMN[CI][2], 0), radius=Radius, fill=False, color='black', linewidth=2.0, \
                           clip_on=False)  # 接收M位置
        plt.gca().add_artist(RXM_C)
        RXN_C = plt.Circle((ABMN[CI][3], 0), radius=Radius, fill=False, color='black', linewidth=2.0, \
                           clip_on=False)  # 接收N位置
        plt1.gca().add_artist(RXN_C)
        # 电阻率画图点用实心圆表示
        Rho_C = plt.Circle((XE, ZE), radius=Radius, fill=True, color='blue', linewidth=2.0, \
                           clip_on=False)  # 拟断面图画图点
    else:
        ME = ABMN[CI, 1]
        NE = ABMN[CI, 2]
        BE = ABMN[CI, 3]
        XE = (NE - ME) / 2 + ME
        ZE = -(BE - NE)
        E_C = plt.Circle((XE, ZE), radius=Radius, fill=False, color='black', linewidth=1.0, \
                         clip_on=False)  # 画图点位置
        plt1.gca().add_artist(E_C)
        # 画出本次测量的地下记录的及发射接收位置
        TXA_C = plt.Circle((ABMN[CI][0], 0), radius=Radius, fill=False, color='red', linewidth=2.0, \
                           clip_on=False)  # 发射A位置
        plt.gca().add_artist(TXA_C)
        TXB_C = plt.Circle((ABMN[CI][3], 0), radius=Radius, fill=False, color='red', linewidth=2.0, \
                           clip_on=False)  # 发射B位置
        plt.gca().add_artist(TXB_C)
        RXM_C = plt.Circle((ABMN[CI][1], 0), radius=Radius, fill=False, color='black', linewidth=2.0, \
                           clip_on=False)  # 接收M位置
        plt.gca().add_artist(RXM_C)
        RXN_C = plt.Circle((ABMN[CI][2], 0), radius=Radius, fill=False, color='black', linewidth=2.0, \
                           clip_on=False)  # 接收N位置
        plt1.gca().add_artist(RXN_C)
        # 电阻率画图点用实心圆表示
        Rho_C = plt.Circle((XE, ZE), radius=Radius, fill=True, color='blue', linewidth=2.0, \
                           clip_on=False)  # 拟断面图画图点

    for CI1 in range(NS):
        if model == 0 or model == 1:
            BE = int(ABMN[CI1, 1])
            ME = int(ABMN[CI1, 2])
            XE = int(ME - BE) / 2 + BE
            ZE = -int(ME - BE)

            E_C = plt.Circle((XE, ZE), radius=Radius, fill=False, color='black', linewidth=1.0, \
                             clip_on=False)  # 画图点位置
        else:
            ME = ABMN[CI1, 1]
            NE = ABMN[CI1, 2]
            BE = ABMN[CI1, 3]
            XE = (NE - ME) / 2 + ME
            ZE = -(BE - NE)
            E_C = plt.Circle((XE, ZE), radius=Radius, fill=False, color='black', linewidth=1.0, \
                             clip_on=False)  # 画图点位置
        #print("E_C == "+str(E_C))
        # print("ABMN[0:0]： ",ABMN[0][0])
        plt.gca().add_artist(E_C)
    plt1.gca().add_artist(Rho_C)
    #plt1.gca().invert_yaxis()
    plt1.gca().set_xlabel('Array Location(m)', size=14)
    plt1.gca().set_ylabel('n-spacing', size=14)
    plt1.gca().xaxis.set_label_position('top')
    plt1.gca().xaxis.tick_top()
    plt1.gca().set_aspect('equal')
    if model == 0 or model == 1:
        abs_yticks=np.arange(-(A_indx_sum_ouji-1), 0)
        print("abs_yticks 1111 = 0/1")
        print(abs_yticks[::-1])
        plt.yticks(abs_yticks[::-1])
    else:
        abs_yticks=np.arange(-(A_indx_sum-1), 0)
        #print("abs_yticks = 3")
        #print(abs_yticks)
        plt.yticks(abs_yticks[::-1])
    # plt1.yticks(np.insert(np.arange(Min_Z, Max_Z + Aspacing / 2, Aspacing / 2), 0, 0))
    plt1.xticks(np.arange(Min_X, Max_X, 1))
    # 用隔离系数代替纵向坐标的标签
    # labels = [item.get_text() for item in ax1.get_yticklabels()]
    # NLAB = len(labels)
    # if model == 0 or model == 1:
    #     labs = np.arange(0, A_indx_sum_ouji, 1)
    # else:
    #     labs = np.arange(0, A_indx_sum, 1)

    # abs_labs = [abs(num) for num in labs]
    # print(abs_labs)
    ax1.set_yticklabels(abs_yticks[::-1])
    # 隐藏图像的底部边界
    ax1.spines['bottom'].set_visible(False)
    fpl = io.BytesIO()

    plt1.savefig(fpl, format="png")
    fio = fpl.getvalue()
    #print("___________________ERTPlot__________________22222")
    plt1.close()

    fpl.close()
    return fio

#频谱
def draw_fft(data,SF,NFFT,filter):
    #print("-------SF,NFFT--------")
    #print(SF,NFFT)

    # SR = int(NFFT) #发射频率
    # #NFFT = 512 #source的长度
    # Phs = 0                             #phase
    # NFFT_calc = int(NFFT/2) + 1                      #number of freqeuncies 频率数
    # #f = SR/2*np.linspace(0,1,NFFT_calc)    #frequencies 频率
    # Delta_T = 1/SR                      #sampling interval
    # A = 1                               #mag
    # t = np.arange(0, NFFT, 1)*Delta_T    #Time vector
    # source = A * signal.square(2*np.pi*SF*t + Phs, duty=0.5)      # square wave
    # FFT =  2 * np.fft.fft(source, n=NFFT)/NFFT
    # f = SR/2*np.linspace(0,1,NFFT_calc)    #frequencies
    SR = NFFT
    if SF == '8' or SF == '4' or SF == '2'or SF == '1':
        SF = int(SF)
    elif SF == '1/2':
        SF = 0.5
    elif SF == '1/4':
        SF = 0.25
    elif SF == '1/8':
        SF = 0.125

    print("====================== SF ======================="+str(SF))
    print("======================filter======================="+str(filter))
    source = []
    filtedData = []
    fs = 1000  # 采样频率
    list_float = data.strip(',').split(",")
    for i in range(len(list_float)):
        # data2 = Decimal(list_float[i])
        data2 = list_float[i]
        #print(data2)
        source.append(float(data2))
    if filter == 1:#低频滤波
        # b, a  =   signal.butter( 50 ,  50 ,  'lowpass' )    #配置滤波器 8 表示滤波器的阶数
        # filtedData  =   signal.filtfilt(b, a, source)   #source为要过滤的信号
        # 设计陷波滤波器来去除50 Hz干扰
        f0 = 50.0  # 噪声的中心频率
        Q = 30.0  # 带宽
        b, a = signal.iirnotch(f0, Q, fs)
        filtedData = signal.filtfilt(b, a, source)
    elif filter == 2:
        b, a  =   signal.butter( 8 , [0.1,0.12],  'bandstop' )    #配置滤波器 8 表示滤波器的阶数
        filtedData  =   signal.filtfilt(b, a, source)   #source为要过滤的信号
    else:
        filtedData = source

    #FFT转化
    FFT =  2 * np.fft.fft(filtedData, n=NFFT)/NFFT
    NFFT_calc = int(NFFT/2) + 1
    f = SR/2*np.linspace(0,1,NFFT_calc)    #frequencies

    plt2.figure(figsize=(15, 7.5))
    plt2.loglog(f,np.abs(FFT[0:NFFT_calc]), 'red')
    plt2.title('Simple Scatterplot')
    plt2.ylabel('Amplitude')
    plt2.xlabel('Frequency(Hz)')
    plt2.ylim([0.00001,10])
    #indx = np.where(f == SF)
    # print('Mag at fundemental is', np.abs(FFT[indx])[0])
    fio = io.BytesIO()

    plt2.savefig(fio, format="png")
    fio1 = fio.getvalue()
    plt2.close()

    fio.close()
    return fio1
#方波
def draw_fft_square_wave(data,SF,NFFT,filter):
    # s = data
    # s = s[1:-1]
    # print(s)
    # source = [float(n) for n in s.split(',')]
    # print(source)
    # s = data
    # s = s[1:-1]
    # source = [Decimal(n) for n in s.split(',')] #str转Decimal

    SR = NFFT
    if SF == '8' or SF == '4' or SF == '2'or SF == '1':
        SF = int(SF)
    elif SF == '1/2':
        SF = 0.5
    elif SF == '1/4':
        SF = 0.25
    elif SF == '1/8':
        SF = 0.125

    print("====================== SF000000000 ======================="+str(SF))
    print("====================== NFFT ======================="+str(NFFT))
    print("======================filter======================="+str(filter))

    source = []
    filtedData = []
    fs = 1000  # 采样频率
    list_float = data.strip(',').split(",")
    for i in range(len(list_float)):
        # data2 = Decimal(list_float[i])
        data2 = list_float[i]
        #print(data2)
        source.append(float(data2))
    if filter == 1:#低频滤波
        # b, a  =   signal.butter( 50 ,  50 ,  'lowpass' )    #配置滤波器 8 表示滤波器的阶数
        # filtedData  =   signal.filtfilt(b, a, source)   #source为要过滤的信号
        # 设计陷波滤波器来去除50 Hz干扰
        f0 = 50.0  # 噪声的中心频率
        Q = 30.0  # 带宽
        b, a = signal.iirnotch(f0, Q, fs)
        filtedData = signal.filtfilt(b, a, source)
    elif filter == 2:
        b, a  =   signal.butter( 8 , [0.1,0.12],  'bandstop' )    #配置滤波器 8 表示滤波器的阶数
        filtedData  =   signal.filtfilt(b, a, source)   #source为要过滤的信号
    else:
        filtedData = source

    #FFT转化
    # FFT =  2 * np.fft.fft(filtedData, n=NFFT)/NFFT
    # NFFT_calc = int(NFFT/2) + 1
    # f = SR/2*np.linspace(0,1,NFFT_calc)    #frequencies
    #plt3.plot(f,np.abs(FFT[0:NFFT_calc]), 'green')#FFT处理

    #print(filtedData)
    Delta_T = 1/int(SR)                     #sampling interval 采样间隔
    t = np.arange(0, NFFT, 1)*Delta_T    #Time vector 时间向量
    plt3.figure(figsize=(15, 7.5))
    plt3.plot(t, filtedData, 'green')
    plt3.ylabel('Amplitude')
    plt3.xlabel('Time (s)')
    # plt3.ylim(-0.1, 0.1)
    plt3.ioff()
    fsq = io.BytesIO()

    plt3.savefig(fsq, format="png")
    fsq1 = fsq.getvalue()
    plt3.close()
    fsq.close()
    return fsq1

def pb2jb(byte_arr):
    """
    python字节码转java字节码
    :param byte_arr:
    :return:
    """
    return [int(i) - 256 if int(i) > 127 else int(i) for i in byte_arr]

def image_to_fig_ax(image_path):
    """
    将图片文件转换为matplotlib的Figure和Axes对象
    参数:
        image_path: 图片文件路径
    返回:
        fig: matplotlib的Figure对象
        ax: matplotlib的Axes对象
    """
    # 读取图片
    img = mpimg.imread(image_path)
    fig = plt.gcf()
    ax = plt.gca()
    # 调整大小（宽度10，高度7）
    fig.set_size_inches(15, 7.5)
    # 关闭自动缩放
    ax.autoscale(False)
    # 创建fig和ax对象
    # fig, ax = plt.subplots(figsize=(15, 7.5))

    # 在ax上显示图片
    ax.imshow(img)

    # 隐藏坐标轴（可选）
    # ax.axis('off')

    return fig, ax
